function [SDX,flag_ok,V]=chol3SDX(parv,flag_check)
%=========================================================================
% [SDX,flag_ok,V]=chol3SDX(parv,flag_check)
%
% Given the 6 elements of the [3 3] Covariance Matrix, V, construct SDX=chol(V) such that
% V=SDX'*SDX while allowing for correlated elements. 
% 
%% *Inputs*
% PARV [6 1] with entried  
% 1) sigma_11 
% 2) sigma_22
% 3) sigma_33
% 4) rho_12
% 5) rho_13
% 6) rho_23 
% 
% FLAG_CHECK   Optional, if ==1 check resulting V matrix is PD
%              default is zero, in which case flag_ok==1 does not mean
%              matrix has been checked 
% 
%% *Outputs*
% SDX           V=SDX'SDX 
% flag_ok       if flag_check=1, means V is P.D. 
% V 
% 
% Alejandro Justiniano Feb 3rd 2011 
%=========================================================================
if nargin < 2 || flag_check~=1; 
    flag_check=0; 
end 
%% Check dimensions and parameters in bounds
flag_ok=0; SDX=[]; 
if length(parv)~=6;error('PARV must be 6x1');end 
if any( parv(1:3) ) < eps; error('First 3 entries of PARV must be standard deviations');end 
if any( abs(parv(4:6)) > 1 ); error('Last 3 entries of PARV must be correlations');end 

%% Unpack 
sigma_11 = parv(1); 
sigma_22 = parv(2); 
sigma_33 = parv(3); 
rho_12   = parv(4); 
rho_13   = parv(5); 
rho_23   = parv(6); 

%% Solve coefficients SDX=[a b c;0 d e;0 0 f]
a=sigma_11; 
b=rho_12*sigma_22; 
c=rho_13*sigma_33; 
d=sigma_22*sqrt( 1-rho_12*rho_12 ); 
e=( rho_23*sigma_22*sigma_33 - b*c )/d; 
ftemp=sigma_33*sigma_33 - e*e - c*c; 
if ftemp < eps; return; end 
f=sqrt(ftemp); 

SDX=zeros(3,3); 
SDX(1,1)=a; 
SDX(1,2)=b; 
SDX(1,3)=c; 
SDX(2,2)=d; 
SDX(2,3)=e; 
SDX(3,3)=f; 
 
if flag_check==0; 
    flag_ok=1; V=[]; 
    return 
end 

%% If flag_check==1, check matrix is PD 
V=(SDX'*SDX);
V=0.5*(V+V'); 
[~,flag_nopd]=chol(V); 
if flag_nopd==0;flag_ok=1;end; 